Gradient Phyto Provinces

Code 
library(dplyr) ; library(tidyr) ; library(ggplot2) ; library(ggpubr) ; library(purrr)

### -- 
generating=F
## -- 

root <- rprojroot::has_file(".git/index")
datadir = root$find_file("data")
funsdir = root$find_file("functions")
savingdir = root$find_file("saved_files")

files_vec <- list.files(funsdir)

for( i in 1:length(files_vec)){
  source(root$find_file(paste0(funsdir,'/',files_vec[i])))
}

###############################
## Loading
###############################
df_geo_abiotics = readRDS(file = paste0(savingdir,'/','df_geo_abiotics'))

dfGrump_longer_filtered = data.table::fread(paste0(datadir,'/','phyto_gradients.csv'))

#list_AitDist = readRDS(paste0(savingdir,'/','list_AitDist'))
list_normalized_geo_abiotics_dists = readRDS(file = paste0(savingdir,'/','list_abio_bio_geo_dist_normallized'))
list_geo_abiotics_dists = readRDS(file = paste0(savingdir,'/','list_abio_bio_geo_dist'))
grid_base = readRDS(file = paste0(savingdir,'/','grid_base'))
Code 
df_evaluation <- readRDS(paste0(savingdir,'/','df_evaluation_choosing_k'))

df_summary = bind_rows(
  df_evaluation %>% 
    group_by(ncluster,method,alpha,DistMetric) %>% 
    summarise_all(.funs = mean) %>% mutate(summary_metric='Mean') %>% ungroup() %>% 
    pivot_longer(cols = -c(ncluster,method,alpha,DistMetric,summary_metric),names_to = 'Metric'),
  
  df_evaluation %>% 
    group_by(ncluster,method,alpha,DistMetric) %>% 
    summarise_all(.funs = sd) %>% mutate(summary_metric='SD') %>% ungroup() %>% 
    pivot_longer(cols = -c(ncluster,method,alpha,DistMetric,summary_metric),names_to = 'Metric')
)

Summary - Biotic - Smetric

Code 
df_summary %>%
  filter(Metric=='avg_within_between_dist2medoid',
         DistMetric=='Bio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Biotic - Sratio')
df_summary %>%
  filter(Metric=='avg_dist_within_dist2medoid',
         DistMetric=='Bio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Biotic - Numerator')
df_summary %>%
  filter(Metric=='avg_dist_between_medoids',
         DistMetric=='Bio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Biotic - Denominator')

Code 
df_summary %>%
  filter(Metric=='avg_within_between_dist2medoid',
         DistMetric=='Geo') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Geo - Sratio')
df_summary %>%
  filter(Metric=='avg_dist_within_dist2medoid',
         DistMetric=='Geo') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Geo - Numerator')
df_summary %>%
  filter(Metric=='avg_dist_between_medoids',
         DistMetric=='Geo') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Geo - Denominator')

Code 
df_summary %>%
  filter(Metric=='avg_within_between_dist2medoid',
         DistMetric=='Abio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Abio - Sratio')
df_summary %>%
  filter(Metric=='avg_dist_within_dist2medoid',
         DistMetric=='Abio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Abio - Numerator')
df_summary %>%
  filter(Metric=='avg_dist_between_medoids',
         DistMetric=='Abio') %>%
  select(-Metric) %>% 
  pivot_wider(id_cols = ncluster:DistMetric,
              names_from = summary_metric ) %>%
  mutate(alpha=factor(alpha)) %>% 
  ggplot(aes(x=ncluster,y=Mean,color=alpha,fill=alpha,linetype=alpha))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  #facet_wrap(alpha~.,nrow = 2)+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Abio - Denominator')

Code 
df_evaluation <- readRDS(paste0(savingdir,'/','df_evaluation_choosing_alpha'))

df_summary = bind_rows(
  df_evaluation %>% 
    group_by(method,alpha,DistMetric) %>% 
    summarise_all(.funs = mean) %>% mutate(summary_metric='Mean') %>% ungroup() %>% 
    pivot_longer(cols = -c(method,alpha,DistMetric,summary_metric),names_to = 'Metric'),

  df_evaluation %>% 
    group_by(method,alpha,DistMetric) %>% 
    summarise_all(.funs = sd) %>% mutate(summary_metric='SD') %>% ungroup() %>% 
    pivot_longer(cols = -c(method,alpha,DistMetric,summary_metric),names_to = 'Metric')
)
Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_within_between_dist2medoid') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Sratio')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_dist_within_dist2medoid') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('S-numerator')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_dist_between_medoids') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('S-numerator')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
df_evaluation <- readRDS(paste0(savingdir,'/','df_evaluation_choosing_alpha_alt'))

df_summary = bind_rows(
  df_evaluation %>% 
    group_by(method,alpha,DistMetric) %>% 
    summarise_all(.funs = mean) %>% mutate(summary_metric='Mean') %>% ungroup() %>% 
    pivot_longer(cols = -c(method,alpha,DistMetric,summary_metric),names_to = 'Metric'),

  df_evaluation %>% 
    group_by(method,alpha,DistMetric) %>% 
    summarise_all(.funs = sd) %>% mutate(summary_metric='SD') %>% ungroup() %>% 
    pivot_longer(cols = -c(method,alpha,DistMetric,summary_metric),names_to = 'Metric')
)
Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_within_between_dist2medoid') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('Sratio')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_dist_within_dist2medoid') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('S-numerator')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
df_summary %>%
  pivot_wider(id_cols = method:Metric,names_from = summary_metric ) %>%
  filter(Metric=='avg_dist_between_medoids') %>% 
  ggplot(aes(x=alpha,y=Mean,color=method,fill=method))+
  geom_line()+
  geom_ribbon(aes(ymin=Mean-SD,ymax=Mean+SD),alpha=0.25)+
  facet_wrap(DistMetric~.,scales = 'free')+
  theme_minimal(base_size = 16)+
  theme(legend.position = 'bottom')+
  ggtitle('S-numerator')+
  scale_x_continuous(breaks = seq(0,1,0.1))

Code 
D = list_normalized_geo_abiotics_dists$bioticDist
vet_alpha = c(0,0.05,0.1,0.2)
D1 = as.dist(D)
D2 = as.dist( 
  (1-as.numeric(vet_alpha[1]))*D +
    as.numeric(vet_alpha[1])*list_normalized_geo_abiotics_dists$geoDist)
D3 = as.dist( 
  (1-as.numeric(vet_alpha[2]))*D +
    as.numeric(vet_alpha[2])*list_normalized_geo_abiotics_dists$geoDist)
D4 = as.dist( 
  (1-as.numeric(vet_alpha[3]))*D +
    as.numeric(vet_alpha[3])*list_normalized_geo_abiotics_dists$geoDist)

####################################   Ward    ####################################
#################################### Creating hclust objs
hclust_0 = hclust(d=D1,method='ward.D')
hclust_0.05 = hclust(d=D2,method='ward.D')
hclust_0.10 = hclust(d=D3,method='ward.D')
hclust_0.20 = hclust(d=D4,method='ward.D')

par(mfrow=c(1,4))
plot(hclust_0   ,labels = F,xlab = '',sub = '')
plot(hclust_0.05,labels = F,xlab = '',sub = '')
plot(hclust_0.10,labels = F,xlab = '',sub = '')
plot(hclust_0.20,labels = F,xlab = '',sub = '')

Code 
par(mfrow=c(1,1))
Code 
nclusters = 7
#################################### With geodist 
mat_cluster_membership_label <- cbind.data.frame(

  ####################################   pam    ####################################
 #pam_a0=cluster::pam(x=D1,k = nclusters,cluster.only = T),
 #pam_a0.05=cluster::pam(x=D2,k = nclusters,cluster.only = T),
 #pam_a0.10=cluster::pam(x=D3,k = nclusters,cluster.only = T),
 #pam_a0.20=cluster::pam(x=D4,k = nclusters,cluster.only = T),
  hclust_a0=cutree(hclust_0   ,k = nclusters),
  hclust_a0.05=cutree(hclust_0.05,k = nclusters),
  hclust_a0.10=cutree(hclust_0.10,k = nclusters),
  hclust_a0.20=cutree(hclust_0.20,k = nclusters)
)
Code 
list_df4plot = df_geo_abiotics %>% bind_cols(
  mat_cluster_membership_label
) %>% 
  tidyr::pivot_longer(cols = c(
    #pam_a0,pam_a0.05,pam_a0.10,pam_a0.20,
    hclust_a0,hclust_a0.05,hclust_a0.10,hclust_a0.20)
    ) %>% 
  mutate(Cluster=as.factor(value),
         method=factor(name,levels=c(
           #'pam_a0','pam_a0.05','pam_a0.10','pam_a0.20',
           'hclust_a0','hclust_a0.05','hclust_a0.10','hclust_a0.20'))) %>%
  group_split(method) #%>% 


list_df4plot %>% purrr::map(
    ~ggplot(., aes(x=Latitude,y=Depth, label = Cluster,fill=Cluster)) + 
      geom_label(size = 3,alpha=0.5) +
      theme_minimal()+
      scale_y_reverse() +
      # scale_colour_gradient2(
      #   low = "deepskyblue", 
      #   mid = "gray", 
      #   high = "orange",
      #   na.value = 'red',
      #   midpoint = median(.$value,na.rm=T)
      # )+
      ggtitle(unique(.$method))+
      theme(legend.position = '')
    #facet_wrap(~ AbioticFactor, labeller = function(x) label_value(x, multi_line = FALSE))
  ) %>% 
  cowplot::plot_grid(plotlist = ., ncol = 4)

Code 
generating=F
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_10neigh_k5'))
  
  clust_member_limits =unlist_coloring_obj(list_cluster_membership_and_bounderies)
  
  plt2 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  saveRDS(plt2,paste0(savingdir,'/','plt2'))
  rm(list_cluster_membership_and_bounderies)
  rm(clust_member_limits)
}
plt2 = readRDS(paste0(savingdir,'/','plt2'))
plt2$limits_faceted

Code 
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_05neigh_k5'))
  
  clust_member_limits = unlist_coloring_obj(list_cluster_membership_and_bounderies)
  rm(list_cluster_membership_and_bounderies)
  
  plt2_n5 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  rm(clust_member_limits)
  saveRDS(plt2_n5,paste0(savingdir,'/','plt2_n5'))
}
plt2_n5 = readRDS(paste0(savingdir,'/','plt2_n5'))
plt2_n5$limits_faceted

Code 
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_mirroredLat_n3_k5'))
  
  clust_member_limits = unlist_coloring_obj(list_cluster_membership_and_bounderies)
  rm(list_cluster_membership_and_bounderies)
  
  plt2_n3 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  rm(clust_member_limits)
  
  saveRDS(plt2_n3,paste0(savingdir,'/','plt2_n3'))
}
plt2_n3 = readRDS(paste0(savingdir,'/','plt2_n3'))
plt2_n3$limits_faceted

Phyto Provinces

Code 
plt2$clustRegion_facetd

Code 
plt2_n5$clustRegion_facetd

Code 
plt2_n3$clustRegion_facetd

Clusters Before matching

Code 
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_10neigh'))
  
  clust_member_limits = unlist_coloring_obj(list_cluster_membership_and_bounderies)
  rm(list_cluster_membership_and_bounderies)
  plt2 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  rm(clust_member_limits)
  saveRDS(plt2,paste0(savingdir,'/','plt2_k10'))
}
plt2 = readRDS(paste0(savingdir,'/','plt2_k10'))
plt2$limits_faceted

Code 
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_05neigh'))
  
  clust_member_limits = unlist_coloring_obj(list_cluster_membership_and_bounderies)
  rm(list_cluster_membership_and_bounderies)
  
  plt2_n5 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  rm(clust_member_limits)
  
  saveRDS(plt2_n5,paste0(savingdir,'/','plt2_n5_k10'))
}
plt2_n5 = readRDS(paste0(savingdir,'/','plt2_n5_k10'))
plt2_n5$limits_faceted

Code 
if(generating){
  list_cluster_membership_and_bounderies = readRDS(
    file =  paste0(savingdir,'/','list_cluster_membership_and_bounderies_mirroredLat_n3'))
  
  clust_member_limits = unlist_coloring_obj(list_cluster_membership_and_bounderies)
  rm(list_cluster_membership_and_bounderies)
  
  plt2_n3 = gen_plots(gridBase = grid_base,clustMemLim = clust_member_limits,df_GeoAbio = df_geo_abiotics)
  rm(clust_member_limits)
  saveRDS(plt2_n3,paste0(savingdir,'/','plt2_n3_k10'))
}
plt2_n3 = readRDS(paste0(savingdir,'/','plt2_n3_k10'))
plt2_n3$limits_faceted

Here I chose \(\alpha=0.01\) and ward.D

Code 
nclusters = 5

mat_cluster_membership_label <- cbind.data.frame(

  ####################################   pam    ####################################
 #pam_a0=cluster::pam(x=D1,k = nclusters,cluster.only = T),
 #pam_a0.05=cluster::pam(x=D2,k = nclusters,cluster.only = T),
 #pam_a0.10=cluster::pam(x=D3,k = nclusters,cluster.only = T),
 #pam_a0.20=cluster::pam(x=D4,k = nclusters,cluster.only = T),
  hclust_a0=cutree(hclust_0   ,k = nclusters),
  hclust_a0.05=cutree(hclust_0.05,k = nclusters),
  hclust_a0.10=cutree(hclust_0.10,k = nclusters),
  hclust_a0.20=cutree(hclust_0.20,k = nclusters)
)


dfPlots = df_geo_abiotics %>%
  bind_cols(mat_cluster_membership_label)

dfGrump_longer_filtered_clustered =  dfGrump_longer_filtered %>% left_join(
  dfPlots %>% transmute(SampleID=SampleID,Cluster = hclust_a0.10)
)

set.seed(1234)
objd_dimRed = distMatrix_dimReduction(
  distMatrix = list_geo_abiotics_dists$bioticDist,
  neighUmap = 5,
  perplexityTsne = 15)

dimReduc_obj <- plot_dimReduc_coords(
  coord_plots = objd_dimRed,
  clusters = dfPlots$hclust_a0.10)

dimReduc_obj$MDS_2d
dimReduc_obj$nMDS_2d
dimReduc_obj$TSNE_2d
dimReduc_obj$umap_2d

Code 
eco_groups_zoop_tax_lat_depth = dfGrump_longer_filtered %>% 
  mutate(Eco_relevant_plank_groups=ifelse(Eco_relevant_plank_groups=='','NotAssigned',Eco_relevant_plank_groups),
         Raw.Sequence.Counts=Corrected_sequence_counts) %>% 
  group_by(SampleID,Eco_relevant_plank_groups) %>% 
  summarise(Raw.Sequence.Counts=sum(Raw.Sequence.Counts)) %>% 
  select(SampleID,Eco_relevant_plank_groups,Raw.Sequence.Counts) %>%
  distinct() %>% 
  group_by(SampleID,Eco_relevant_plank_groups) %>% 
  summarise(Sum_RawCounts = sum(Raw.Sequence.Counts)) %>% ungroup %>% 
  tidyr::pivot_wider(
    id_cols = SampleID,
    names_from = Eco_relevant_plank_groups ,
    values_from = Sum_RawCounts,
    values_fill = 0.01) %>%
  data.frame() %>% mutate(SampleID=factor(SampleID)) %>% 
  mutate(across(where(is.numeric))/rowSums(across(where(is.numeric)))) %>% 
  arrange(SampleID)

names_eco_relevant = colnames(eco_groups_zoop_tax_lat_depth[,-1])

eco_relevant_plot <- 
  dfPlots %>% mutate(SampleID=factor(SampleID)) %>% 
  left_join(eco_groups_zoop_tax_lat_depth) 
  
eco_relevant_plot_df = 
  eco_relevant_plot %>%
  transmute(
    SampleID,Latitude,Depth,Cluster=hclust_a0.10
  ) %>% 
  left_join(
    eco_relevant_plot %>% 
      select(SampleID, one_of(names_eco_relevant))
  )
    

eco_relevant_plot_df %>%
    tidyr::pivot_longer(-c('SampleID','Latitude','Depth','Cluster'),
                        names_to = 'Eco_Relevant',values_to = 'RA') %>% 
  mutate(Cluster=factor(Cluster)) %>% 
    ggplot2::ggplot(aes(x=Latitude,y=Depth,size=RA,color=Cluster)) +
    geom_point(alpha=0.5)+
    facet_wrap(Eco_Relevant~.)+
    scale_y_reverse()+
    theme_minimal(base_size = 16)+
    theme(legend.position = 'bottom')